The COLA Anomaly Coupled Model Global SST Predictions

contributed by Ben P. Kirtman^1,2 and Dughong Min²

¹George Mason University and ²Center for Ocean-Land-Atmosphere Studies

The anomaly coupled GCM (ACGCM) combines the COLA atmospheric GCM and the Geophysical Fluid Dynamics Laboratory (GFDL) Modular Ocean Model (MOM), version 3.0. Brief descriptions of these models and the coupling procedures are given below. Details of how well the model performs in long climate simulations is described in detail in Kirtman et al. (2002) and Kirtman and Shukla (2002). An in depth discussion of the hindcast skill is provided in Kirtman (2002) and Kirtman (2003).

a. Atmospheric Model

The dynamic core used in the National Center for Atmospheric Research (NCAR) Community Climate Model (CCM) version 3.0 has been adopted (Schneider, 2001). The dynamic core is spectral (triangular truncation at total wavenumber 42) with semi-Lagrangian transport. There are 18 unevenly spaced -coordinate vertical levels. The parameterization of the solar radiation is after Briegleb (1992) and terrestrial radiation follows Harshvardhan et al. (1987). The deep convection is an implementation of the Relaxed Arakawa-Schubert (RAS) scheme of Moorthi and Suarez (1992) described by DeWitt (1996). The convective cloud fraction follows the scheme used by the NCAR CCM (Kiehl et al., 1994; see DeWitt and Schneider, 1996 for additional details). The model includes a turbulent closure scheme for the subgrid scale exchange of heat, momentum, and moisture after Miyakoda and Sirutis (1977) and Mellor and Yamada (1982). Additional details regarding the AGCM physics can be found in Kinter et al. (1988) and DeWitt (1996). Model documentation is given in Kinter et al. (1997).

b. Ocean Model

The ocean model is version 3 of the GFDL MOM (Pacanowski and Griffies, 1998), a finite-difference treatment of the primitive equations of motion using the Boussinesq and hydrostatic approximations in spherical coordinates. The domain is that of the world ocean between 74^oS and 65^oN. The coastline and bottom topography are realistic except that ocean depths less than 100 m are set to 100 m and the maximum depth is set to 6000 m. The artificial high-latitude meridional boundaries are impermeable and insulating. The zonal resolution is 1.5^o. The meridional grid spacing is 0.5^o between 10^oS and 10^oN, gradually increasing to 1.5^o at 30^oN and 30^oS and fixed at 1.5^o in the extratropics. There are 25 levels in the vertical with 17 levels in the upper 450 m. The vertical mixing scheme is the non-local K-profile parameterization of Large et al. (1994). The horizontal mixing of tracers and momentum is Laplacian. The momentum mixing uses the space-time dependent scheme of Smagorinsky (1963) and the tracer mixing uses Redi (1982) diffusion along with Gent and McWilliams (1990) quasi-adiabatic stirring.

c. Coupling Strategy

The anomaly coupling strategy is described in detail in Kirtman et al. (1997) and in Kirtman et al. (2002). The main idea is that the ocean and atmosphere exchange predicted anomalies, which are computed relative to their own model climatologies, while the climatology upon which the anomalies are superimposed is specified from observations. The anomaly coupling strategy requires atmospheric model climatologies of momentum, heat and fresh water flux, and an ocean model SST climatology. Similarly, observed climatologies of momentum, heat and fresh water flux and SST are also required. The model climatologies are defined by separate uncoupled extended simulations of the ocean and atmospheric models. In the case of the atmosphere, the model climatology is computed from a 30 year (1961-1990) integration with observed specified SST and sea ice. This SST is also used to define the observed SST climatology. In the case of the ocean model SST climatology, an extended uncoupled ocean model simulation is made using 30 years of 1000 mb National Centers for Environmental Prediction (NCEP) reanalysis winds. As with the SST, this observed wind stress product is used to define the observed momentum flux climatology. The heat flux and the fresh water flux in this ocean-only simulation is parameterized using damping of SST and sea surface salinity to observed conditions with a 100 day time scale. The heat and fresh water flux "observed" climatologies are then calculated from the results of the extended ocean only simulation. The ocean and atmosphere model exchange daily mean fluxes and SST once a day.

d. Retrospective Forecast Experiments

In order to assess the potential predictive skill of the coupled model, a large sample of retrospective forecast experiments have been made and compared to available observations. The retrospective forecasts or hindcasts cover the period 1980-1999. A twelve-month hindcast is initialized each January, April, July and October during this 20-year period. For each initial month, an ensemble of six hindcasts is run, yielding a total of 480 retrospective forecasts to be verified. The ocean and atmosphere initial states and the method of generating the ensemble members are described below.

The ocean initial conditions were taken from a 1980-1999 ocean data assimilation produced at GFDL using variational optimal interpolation (Derber and Rosati, 1989). The GFDL ocean initial states were generated using a somewhat higher resolution ocean model than that used in the ocean component of the ACGCM with identical physics and parameter settings. In the forecast experiments presented here, the ocean initial states were interpolated to the lower resolution.

The atmospheric initial states are taken from an extended atmosphere-only simulation with observed prescribed SST. The atmospheric ensemble members were obtained by resetting the model calendar back one week and integrating the model forward one week with prescribed observed SST. In this way, it is possible to generate an unlimited sample of initial conditions that are synoptically independent (separated by one week) but have the same initial date.

e. Deterministic Skill Assessment

Typically, the most common skill assessment is to calculate the Niño3.4 (or Niño3) anomaly correlation and root mean square error (rmse). These are often referred to as the Niño3.4 skill scores. Figure 1 shows the skill scores for Niño3.4 in a slightly different format. The skill scores for the ensemble mean hindcast are depicted by the filled circles. The shaded region indicates the range (plus or minus one standard deviation) of these skill scores calculated from 10,000 possible one-member ensembles that were randomly chosen. In other words, this is the skill score for 10,000 possible combinations of hindcasts where there is one hindcast for each initial month. This is a small sub-sample when one considers that fact that there are 6⁸⁰ possible one member ensembles. A sub-sample of 5,000 shows similar spread in the correlation, indicating that this is a robust estimate of the uncertainty in the correlation.

The skill scores for a persistence hindcast are also shown in Fig. 1. In terms of the correlation coefficient, both the ensemble mean and the one-member ensembles significantly beat a persistence hindcast for all lead times. The ensemble mean correlation coefficient remains above 0.6 for lead times up to nine months, and, for a one-member ensemble, the correlation coefficient remains above 0.6 for lead times up to 7-9 months. In the first two months, the rmse for persistence, the ensemble mean and the one-member ensembles are indistinguishable from each other. Beyond three months lead time, the CGCM hindcasts (ensemble mean or one-member) beat persistence. At a lead time of seven months, the persistence hindcast rmse becomes larger than the climatological hindcast rmse (0.88) indicating that, for lead times greater than seven months, persistence is a particularly poor measure of minimum skill.

The shaded region in Fig. 1 is plotted for two reasons. First, it provides a conservative estimate of the range of uncertainty in these skill scores. We have not provided error estimates for the ensemble mean skill scores, but suggest that the error bar would be less than that for the one-member ensemble. Second, the shaded region also indicates the increase in skill that is expected by using an ensemble of hindcasts produced by a single model. In this case, the increase in skill comes from averaging across an ensemble of six members, and for both the correlation coefficient and the rmse, the ensemble average is about one standard deviation better than the expected value for a one-member ensemble.

f. Probabilistic Skill Assessment

ROC curves for the Niño3.4 hindcasts at lead times of three, six and twelve months are shown in Fig. 2. In calculating the ROC curves, the corresponding contingency tables have been aggregated over all model grid points in the Niño3.4 region. There is a contingency table for each grid box and the aggregation consists of summing the respective elements of the contingency for all grid boxes in the Niño3.4 region. For all three lead times, there are three curves or equivalently we have considered three different events: (i) warm events (upper tercile), (ii) cold events (lower tercile) and (iii) near normal (middle tercile), where both the retrospective forecasts and the observations have been normalized by their local standard deviation. This ability to easily verify the hindcast skill of warm events and cold events separately is one of the advantages of the ROC calculation.

An ideal probabilistic forecast system would have relatively large hit rates and small false alarm rates so that all the points on the ROC curve would cluster in the upper left corner of the diagram. For a relatively poor forecast system, all the points of the ROC curve would lie very close to the dashed diagonal line indicating that the hit rate and the false alarm rate were nearly the same (i.e., no skill). The six interior points on the ROC curve indicate number of ensemble members forecasting a particular event. Progressing from the lower left corner to the upper right corner, the first interior point on the curve (i.e., the first point does not lie on the lower left corner) corresponds to six out of six ensemble members forecasting the event. This point indicates how skillful the forecast system is when the model consistently forecasts a given event to occur. Not surprisingly, the hit rates are modest, but if the false alarm rate is also low, then a confident forecast is very useful. The second point indicates that five out of six ensemble members forecasted the event and the remaining points along the curve vary similarly. The last interior point corresponds to only one ensemble member forecasting the event and the hit rates tend to be high. Unfortunately, the false alarm rates are also high so that there is considerable risk when taking action based on only one ensemble member forecasting an event.

For lead times of both three and six months (upper and middle panel), both warm and cold events are fairly well predicted. The false alarm rates are low and the hit rates are relatively high when the agreement among the ensemble members is relatively large. The model is not particularly overconfident, and, thus, there are no serious limitations associated with the ensemble strategy. In this case, overconfidence would be apparent by large false alarm rates for high forecast agreement. Even when there is only modest agreement among the ensemble members, the hit rates are significantly larger than the false alarm rates. For a near normal forecast, the three month lead has some skill although smaller than the extemes, whereas for both six and twelve month leads, the ROC curve lies close to the diagonal indicating little skill. At twelve months lead time, there is a considerable drop in skill. High confidence forecasts for warm events are only marginally better than those for near normal, suggesting that a confident forecast for a warm event at twelve months lead time is not particularly useful.

g. Forecasts

Forecast for NINO3.4 SSTA initialized in January 2003 and in April 2003 are presented in Fig. 3. All the forecasts indicate a cooling trend through the end of 2003. The January forecasts (upper panel) indicate a non-zero probability of near normal conditions for late 2003, whereas the April forecasts agree on cold temperatures for September-December 2003. The forecast initialized in April (lower panel) suggest a probable return to near normal in early 2004, but two ensemble member predict that the cold temperatures will persist until February 2004.

The most probable (i.e. ensemble mean) SSTA for the forecast initialized in April 2003 is shown in Fig. 4. The top panel shows the ensemble mean forecast for June-August 2003, the middle panel show the ensemble mean forecast for September-November 2003 and the bottom panel shows the ensemble mean forecast for December 2003 - February 2004. The shading in the figure indicates the percentage of ensemble members predicting either the upper or lower tercile. If the ensemble mean is positive (negative), the shading indicates the percentage of ensemble members predicting the upper (lower) tercile.

The most probable forecast calls for relatively cold conditions in the eastern tropical Pacific to persist and, perhaps even strengthen, through February 2004. The consistency among the ensemble members steadily decrease through the forecast period, so that only 4 out of 6 ensemble members are forecasting cold conditions in the eastern Pacific by DJF2004. Near normal conditions prevail in the tropical Atlantic and Indian Oceans. There are strong ensemble mean SSTA in the north Pacific. These anomalies are consistent among all the ensemble members, and in hindcast mode with this level of consistency the model is skillful.

References:

Barnston, A. G., M. Glantz, Y. He, 1999: Predictive skill of statistical and dynamical climate models in SST forecasts during the 1997-98 El Nino and the 1998 La Nina onset. Bul. Amer. Met. Soc., 217-243.

Briegleb, B. P., 1992: Delta-Eddington approximation for solar radiation in the NCAR community climate model. J. Geophys. Res., 97, 7603-7612.

Derber, J., and A. Rosati, 1989: A global oceanic data assimilation system. J. Phys. Oceanogr., 19, 1333-1347.

DeWitt, D. G., 1996: The effect of the cumulus convection on the climate of the COLA general circulation model. COLA Tech. Rep. 27, 69 pp. [Available from COLA, 4041 Powder Mill Road, Suite 302, Calverton, MD 20705.]

DeWitt, D. G., and E. K. Schneider, 1996: The Earth radiation budget as simulated by the COLA GCM. COLA Tech. Rep. 35, 39 pp. [Available from COLA, 4041 Powder Mill Road, Suite 302, Calverton, MD 20705.]

Gent, P. R., and J. C. McWilliams, 1990: Isopycnal mixing is ocean circulation models. J. of. Phys. Oceanography, 25, 150-155.

Harshvardhan, R. Davis, D. A. Randall and T. G. Corsetti, 1987: A fast radiation parameterization for general circulation models. J. Geophys. Res., 92, 1009-1016.

Kiehl, J. T., J. J. Hack and B. P. Briegleb, 1994: The simulated Earth radiation budget of the National Center for Atmospheric Research community climate model CCM2 and comparisons with the Earth Radiation Budget Experiment (ERBE). J. Geophys. Res., 99, 20815-20827.

Kinter, J. L. III, J. Shukla, L. Marx and E. K. Schneider, 1988: A simulation of winter and summer circulations with the NMC global spectral model. J. Atmos. Sci., 45, 2468-2522.

Kinter, J. L. III, D. G. DeWitt, P. A. Dirmeyer, M. J. Fennessy, B. P. Kirtman, L. Marx, E. K. Schneider, J. Shukla and D. Straus, 1997: The COLA Atmosphere-Biosphere General Circulation Model Volume 1: Formulation. COLA Tech. Rep. 51, 46pp. [Available from COLA 4041 Powder Mill Road, Suite 302, Calverton, MD 20705.]

Kirtman, B. P., 2003: The COLA anomaly coupled model: Ensemble ENSO prediction. Mon. Wea. Rev., (in press).

Kirtman, B. P., 2002: The COLA anomaly coupled model: Ensemble ENSO prediction. COLA Tech. Rep. 125, 52 pp. [Available from COLA 4041 Powder Mill Road, Suite 302, Calverton, MD 20705; http://www.iges.org/pubs/tech.html ].

Kirtman, B. P., Y. Fan and E. K. Schneider, 2002: The COLA global coupled and anomaly coupled ocean-atmosphere GCM. J. Climate, 15, 2301-2320.

Kirtman, B. P., and J. Shukla, 2002: Interactive coupled ensemble: A new coupling strategy for CGCMs. Geophys. Res. Let., 29, 17-20.

Kirtman, B. P., J. Shukla, B. Huang, Z. Zhu and E. K. Schneider, 1997: Multiseasonal prediction with a coupled tropical ocean global atmosphere system. Mon. Wea. Rev., 125, 789-808.

Large, W. G., J. C. McWilliams, and S. C. Doney, 1994: Oceanic vertical mixing: A review and a model with a nonlocal boundary layer parameterization. Rev. of Geophys., 32, 363-403.

Mellor, G. L., and T. Yamada, 1982: Development of a turbulence closure model for geophysical fluid problems. Rev. Geophys. Space. Phys., 20, 851-875.

Miyakoda, K., and J. Sirutis, 1977: Comparative integrations of global spectral models with various parameterized processes of sub-grid scale vertical transports. Beitr. Phys. Atmos., 50, 445-480.

Moorthi, S., and M. J. Suarez, 1992: Relaxed Arakawa-Schubert: A parameterization of moist convection for general circulation models. Mon. Wea. Rev., 120, 978-1002.

Pacanowski, R. C., and S. M. Griffies, 1998: MOM 3.0 Manual, NOAA/Geophysical Fluid Dynamics Laboratory, Princeton, USA 08542.

Redi, M. H., 1982: Oceanic isopychnal mixing by coordinate rotation, J. of Phys. Oceanogr., 12, 1155-1158.

Schneider, E. K., 2001: Understanding differences between the equatorial Pacific as simulated by two coupled GCM's. J. Climate, 15, 449-469.

Smagorinsky, J., 1963: General circulation experiments with the primitive equations: I. The basic experiment. Mon. Wea. Rev.,91, 99-164 91,99-164.

Figure 1: Niño3.4 Skill scores: (a) Anomaly correlation coefficient and (b) root mean square error. The skill of the persistence forecast is denoted by the solid blue line and the ensemble mean skill is indicated by the red line with filled circles. The green shaded region indicates ± for 10,000 randomly chosen 1-member ensembles.

Figure 2: Niño3.4 ROC curves for a lead time of three months (top panel), six months (middle panel) and twelve months (bottom panel). Warm events (upper tercile) are denoted in red, near norm conditions (middle tercile) are denoted in green and cold events (lower tercile) are denoted in blue.

Figure 3: NINO3.4 SSTA evolution for forecast initialized in January 2003 (upper panel) and April 2003 lower panel. The ensemble mean is denoted by the dot-dashed curve. The individual ensemble members are denoted by the solid curves. The available observations are denoted by the thick solid curve with open circles.

Figure 4: Predicted ensemble mean SSTA (contours) using April 2003 initial conditions. The shading corresponds to the percentage of ensemble members predicting in either the upper or lower tercile. If the ensemble mean is positive (negative), the shading indicates the percentage of ensemble members predicting upper (lower) tercile. The contour interval for SSTA is 0.5^oC.